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Abstract 

We propose a random matrix model that interpolates between the chiral random 
matrix ensembles and the chiral Poisson ensemble. By mapping this model on a non- 
interacting Fermi-gas we show that for energy differences less than a critical energy E c 
the spectral correlations are given by chiral Random Matrix Theory whereas for energy 
differences larger than E c the number variance shows a linear dependence on the energy 
difference with a slope that depends on the parameters of the model. If the parameters are 
scaled such that the slope remains fixed in the thermodynamic limit, this model provides 
a description of QCD Dirac spectra in the universality class of critical statistics. In this 
way a good description of QCD Dirac spectra for gauge field configurations given by a 
liquid of instantons is obtained. 

PACS: 11.30.Rd, 12.39.Fe, 12.38.Lg, 7I.30.+h 
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1 Introduction 



By now it has been well-established that the smallest eigenvalues of the QCD Dirac 
operator are correlated according to a Random Matrix Theory with the global symmetries 
of the QCD partition function [0J (see 0, |3| for recent reviews and a complete list of 
references). In particular, this has been confirmed by the analysis of the low-energy 
effective theory JI], ^ ^, [/J, universality studies || [| |ltj, [IT], |I2| , lattice QCD simulations 
lj, |, 0, 0, |TJ, |Cj, ^1] g2[ and the study two-sublattice theories with disorder 



SB], E3L EH EBf. This means that the dynamical details of QCD are not important on 



energy scales of the order of the average level spacing. The natural question that can be 
asked is at what energy scale the dynamics of QCD becomes relevant and how does this 
manifest itself in the Dirac spectrum. 

The answer to this question has been understood within the context of effective theories 
p7| , |28| , [|, |], |6], 0. The effective theory for the QCD Dirac spectrum is known as the 
partially quenched effective partition function and was originally introduced to study the 
quenched approximation in QCD |2"U| , |50| . The central observation is that in the domain 
where the mass dependence of the effective partition function is given by the contribution 
of the constant fields (the zero momentum modes) the Dirac eigenvalues are correlated 
according to chiral Random Matrix Theory. The relevant mass scale can thus be identified 
as the scale for which the Compton wavelength of the lightest particle becomes equal to 
the size of the box i.e., M ~ l/L s . In QCD, in the phase of spontaneously broken chiral 
symmetry, the lightest particles are the Goldstone modes with a mass given by M 2 ~ Km 
(with m the quark mass and K = S/F 2 in terms of the pion decay constant F and the 
chiral condensate £). The critical scale is thus given by 



In the context of disordered condensed matter systems |3TJ, |32], [IjJ this energy scale is 
known as the Thouless energy and also in this article we will adopt this name. 

A more intuitive interpretation of the Thouless energy has been given in the theory of 
mesoscopic systems JJTJ. The time scale h/E c is the time for which an initially localized 
wave packet diffuses all over space. For this reason the eigenvalues are correlated ac- 
cording to Random Matrix Theory for energy differences below E c (known as the ergodic 
regime). At shorter time scales, different wave functions do not necessarily overlap result- 
ing in a weakening of correlations of the corresponding eigenvalues. For energy differences 
beyond the inverse elastic collision time r e the corresponding eigenvalues are completely 
uncorrelated (the Poisson ensemble). The domain inbetween E c and h/r e is known as 
the diffusive or Altshuler-Shklovskii domain. A third energy scale is the average level 
spacing AA. The ratio E C /AX is identified in mesoscopic physics as the dimensionless 



2 



conductance. It is equal to the number of subsequent levels correlated according to Ran- 
dom Matrix Theory. The existence of these domains has been confirmed by numerical 



simulations of the Anderson model [34 



In QCD the average level spacing is related to the order parameter of the chiral 
phase transition, the chiral condensate, by the Banks-Casher formula |35| according to 



AA = it/Ty (with V the volume of Euclidean space-time). The prediction from ([!]) is 
that the number of eigenvalues correlated according to chRMT is of the order F 2 \/V. 

For increasing disorder the number of subsequent eigenvalues described by Random 
Matrix Theory decreases. For strong disorder we expect that all states become localized 
with uncorrelated eigenvalues. We thus expect a critical value of the disorder for which 
the three scales, AA, E c and h/r e coincide. In particular, the dimensionless conductance 
becomes volume independent j3f|. It has been conjectured |J7| that at this point the 



eigenvalue correlations are described by a new universality class known as critical statistics 
for a review). In this class, only the short range correlations of the eigenvalues 



.see 



are described by the usual random matrix ensembles whereas the number variance, S 2 (n), 
shows a linear n-dependence beyond this domain. What is relevant for QCD is that such 
behavior has been observed in numerical simulations of the 4-dimensional Anderson model 
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The volume dependence of the Thouless energy has been investigated by means of lat- 
tice QCD simulations |16|, 0, |22|] and instanton-liquid simulations . In essence, results 



from lattice QCD simulations are in complete agreement with theoretical results from 
partial quenched chiral perturbation theory. However, the results from instanton liquid 
simulations seem to deviate from the prediction ([!]) with a scale independent constant 
K\ in that case the Thouless energy only shows a weak volume dependence. This raises 
the question whether the Dirac eigenvalues might be described by critical statistics. To 
address this issue we generalize a random matrix model for critical statistics [11] to 



include the chiral symmetry of the QCD partition function (section 2). In section 3 we 
map our model on a partition function of noninteracting fermions. This model is solved 
in the semi-classical limit in section 4, where we obtain analytical expressions for the 
microscopic spectral density and the two-point correlation function. Comparisons with 
instanton simulations are shown in section 5 and concluding remarks are made in section 
6. 

2 Definition of the Model 

The random matrix model of Moshe and Neuberger and Shapiro [[5J is defined by the 
partition function 

Z = I dHe-™ H ' J dUe -bTr([u,H][u, H ]f)^ (2) 
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where if is a Hermitian and U a Unitary N xN matrix. The integration measures dH and 
dll are given by the Haar measure. This model can be interpreted as the zero-dimensional 



limit of the Kazakov-Migdal model J42|. It interpolates between the Gaussian Unitary 
Ensemble (b = 0) and the Poisson Ensemble (b — > oo). Using the invariance of the 
measure, the integral over U can be replaced by an integral over the eigenvalues of U. For 
b — ► oo this partition function is dominated by matrices H that commute with arbitrary 
diagonal unitary matrices. This set of matrices is the ensemble of diagonal Hermitian 
matrices which is known as the Poisson ensemble. What is nice about this model is 
that it preserves the unitary invariance which enables us to take full advantage of the 
existing random matrix theory methods. In order to obtain a nontrivial 6-dependence in 
the thermodynamic limit, the parameter b has to be scaled as 

b = h 2 N 2 . (3) 



It has been shown that in this limit the model (H) is equivalent f4l| to both a banded 



random matrix model JEJ with a power-like cutoff |44| and to random matrix models 
with a Tr log 2 H probability potential j|5| . The correlation functions of the latter model 



have been derived by means of g-orthogonal polynomials |46|] and Painleve equations 



57J. 




In this paper we are interested in chiral random matrix ensembles defined as ensembles 
of N x N random matrices with the structure 

(4) 

where C is an arbitrary complex matrix n x (n+u) matrix (N = 2n+u). Since the matrix 
D has exactly v zero eigenvalues v is interpreted as the topological quantum number. The 
chiral Gaussian Unitary Ensemble (chGUE) with Nf flavors is defined as the ensemble 
of matrices D with matrix elements C distributed according to the Gaussian probability 
distribution 

P(C) ~ det"' (D + m)e~ NTrCC ' . (5) 

Here, for simplicity we have taken all quark masses equal to m. The probability distribu- 
tion P(C) has the unitary invariance 

C UCV-\ (6) 

where U and V are unitary matrices. Since an arbitrary complex matrix can always be 
brought to diagonal form by this transformation this invariance allows us to factorize the 
probability distribution in a product over the eigenvalues of C and the unitary matrices 
that diagonalize C. 
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The generalization of the model of Moshe, Neuberger and Shapiro to the chiral en- 
sembles is immediate. The interpolating model is defined by the partition function 

Z V = J dCdet N f(D + m)e-^ TlD ' D J dUe~^^ D ^ D ^\ (7) 

where U has the chiral block structure 



«= n I I (<»> 




with [/aiinxn unitary matrix and V an (n + v) x (n + v) unitary matrix. If E and h 
are constants as in fl7|) we will refer to this model as the critical chiral unitary ensemble. 
As we will see below, in order to make contact with the Thouless energy in the partially 
quenched effective partition function we have to scale h with an additional factor 1/x/iV. 
The unitary invariance of this partition function follows from the invariance of the Haar 
measure. In comparison to | iD| , an additional factor ~ E 2 /h has been included in the 
probability distribution of the matrix elements. As we will see below, this will guarantee 
that in the thermodynamic limit the spectral density is /i-independent to leading order in 
h. We will also find that the partition function is normalized such that the E represents 
the chiral condensate by means of the Banks-Casher relation E = 7rp(0)/N (with p(0) the 
spectral density around A = 0). 

Decomposing into the blocks of D and U the partition function can be written as 

Z u = m» j dCdet N > [C^C + m 2 ]e -(i +2 /^)iW j dUdVe zW^ucv-^ (g) 
The arbitrary complex matrix C can be decomposed according to 

C = UtAU 2 , (10) 

and the integral over C can be expressed as an integral over the eigenvalues A& and the 
unitary matrices JJ-y and U%. Up to a irrelevant constant, the Jacobian of this transfor- 
mation is given by 

j(A)=A 2 ({A 2 })nAr +i - (ii) 

k 

where 

A({A. 2 }) = L[(^-^) (12) 

k<l 

is the Vandermonde determinant. Because of the unitary invariance, the U± and Ui 
dependence in the second exponent can be absorbed in a redefinition of [Zand V, and the 
integrations over U\ and U2 just result in an overall constant. Remarkably, the integral 
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over U and V in (H) is an Itzykson-Zuber type integral which is known analytically 



48, 49, 50 



f JTT^/^^ukv^k _ n detfc,; \I u (z\ k \i)\ 

I dudVe ~ °*mmuw- (13) 

The Vandermonde determinants cancel resulting in the partition function 

/n n 2 

II d\ k II (Al + m 2 ff det (A fc A0^-^ +/l2iv2 ^^ + Af) 4( ^ 2iv2AfcA ^ 
1111 



k=l k=l 



(14) 



The instanton-liquid Dirac spectra that will be described by this model were obtained in 
the quenched approximation and for zero total topological charge. Therefore we will not 
attempt to solve this model for arbitrary Nf but instead focus on the technically simpler 
case of Nf = 0. Since the topological charge does not give rise to additional complications 
we will consider the case of arbitrary v. Below we will thus analyze the joint probability 
distribution 



MAi, • • • , A n ) = det 

k,l 



(A jfc A z )le-(i + ' l2jv2 )ir^+^) 4( ^ 2iv2Afe ^ : 



(15) 



3 Fermi-Gas Representation of the Interpolating chi- 
ral Random Matrix Model 



In this section we rewrite the joint probability distribution (|15f) in terms of the fermionic 
n-particle matrix element 

p u {x!, ■■■,x n ) = C(xi ■ ■ -x n \e~ pH \xi ■ -x n ), (16) 

where H is the separable Hamiltonian 

H = Y.(-dl + ^-^ + u 2 xl) ) (17) 

k 

and C is an irrelevant constant. The eigenfunctions of the single particle Hamiltonian are 
known in terms of Laguerre polynomials. Specifically, 

{-d 2 x + 4 ^ ~ 1 + tu 2 x 2 )(f) n = uj n (f) n , (18) 



has the solutions 



Un = (An + 2u + 2)u, (19) 



r - - (x^) V+1 Ll(x 2 u) . (20) 
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The normalization factor h n = (n + i>)\/n\ ensures that the n are normalized to unity 

i 

4>n{x)<f)n{x)dx = I. (21) 

Using completeness, the many-particle matrix element can be written as 

{x x ■ ■ ■ x n \eT m \x x ■ ■ ■ x n ) = det V ^„(x i )e- / * Wn n (a; j ) (22) 

ij 

J n 

The sum over n can be performed analytically using the identity 

^ L v n {x)L" n {y)z n _ {xyz)-»l 2 / _ x + y\ ( 2^xyz\ 
to K ~ l-z 6XP l Z l-z) lv \l-z)- 

With z = exp(— 4/9o;) this results in 



(23) 



I *^ i * * * *^ tt, ^ — d.G"fc 

ij 



ui^/XiXj ( uj cosh 2(5uj . 2 2 \ ( u 



sinh 2{3u 6XP ^ 2 sinh 2f3u [%l + Xj) ) v ^sinh 2/3u> XiXj , 
Comparing this expression with fll5|) we find that 



(24) 



u = V4h 2 N 2 + 1 — , (25) 
2h 

cosh 2^ = 1+2^- (26) 

The joint probability distribution of the eigenvalues is thus given by an n— particle 
diagonal matrix element of the density operator. The average spectral density is equal to 
the average particle density. It is obtained by integrating over the positions of all particles 
except one. The integral can be performed by rewriting the matrix elements (|22|) in terms 
of a sum over permutations ir and 7r', 

— E E CT W CT (^0^ (1) (zi)---^ (n) W (27) 

where cr(-) is the sign of the permutation. Performing the integrations over X2, ■ ■ ■ , x n , 
by orthogonality we find that the only nonzero contribution is for 7r = ir'. Then the 
remaining sum over ir(ui) is just a sum over uii, ■ ■ ■ ,u n . For the canonical ensemble we 
thus find the one-particle density, 

<A*)) = (E*(*-*k)) = i- E EI^(^)l 2e ^ Efe ^- (28) 

k n uji < ■ ■ ■ <wjv i 

In an occupation number representation this can be rewritten as 

</>(*)> = E E^toiv^^, (29) 
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where the sum is over all rii G {0, 1} subject to the condition given below the summation 
sign and the canonical partition function is given by 

Z n = e~^ n ^\ (30) 

For completeness we give the following exact expressions for the canonical partition 
function fH], |52j 

e -2(3u)n 2 -2f3uu)n 

Zn = (1 - e" 4 ^)(l - e- 8 ^) ■■■{!- e~^y ^ 



and the one-particle density 



N 

(p(x)) = Y,\M^)\ 2 Y,(- 1 Y e '^ w{kr ~ nr+1 ' r{r+l)) II (l-e~ 4/3 " p ). (32) 

k r=0 p=N~r 

They have been obtained by writing the constraint Y^i n i — n as 

5(N-J2n*) = ^-f d9e iQ{ - N -^i ni \ (33) 

ZTT J 



and summing the geometric series after using the explicit expression ( |20| ) for the io n . 

Instead of working with the canonical ensemble we eliminate the constraint on the sum 
over the rij by working with the grand canonical ensemble. The grand canonical partition 
function is defined by 

Z = J2z n Z n = H(l + e-^+^), (34) 

n k 

where z = is the fugacity. The one particle density in the grand canonical ensemble 
is given by 

(p(x)} = Wz n Z n p(x) 

L n 

= j E Y, n i\H^)\ 2 e-^ nk ^\ (35) 

ni£{0,l} i 

which can be evaluated to be 

The fugacity z = e 13 ^ is determined by the condition that the total number of particles is 
equal to n, i.e. 

E ~, \, V = n. (37) 
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The connected two-particle correlation function is defined by 



Ra(x,y) = ^2(S(x - x k )5{y - y^) - (p(x))(p(y)) 



(38) 



It can be obtained from the n-particle matrix element (^) by integration over all coor- 
dinates (eigenvalues) except x\ and x 2 . The remaining sum over u^i) and coV( 2 ) can be 
rewritten as 



E E- 



[X 2 ~ (puAXl)<pwA x 2> 



0J k 



(39) 



The diagonal term, i = j, is zero allowing us to include it in the summation. The 
first term can then be identified as the square of the average particle density. It is the 
disconnected contribution to the two-particle correlation function. The connected part of 
the two-particle distribution function can thus be written as 



R 2 (x 1 ,x 2 ) = -— E<M x i> 

u>i<—<Wn i,j 



(40) 



This correlation function does not include the contributions from the self-correlations 
of the eigenvalues. After all, our starting point was the joint distribution of different 
eigenvalues. In an occupation number representation R 2 (xi,x 2 ) simplifies to 



R 2 (x 1 ,x 2 ) 



■T £ 



E 



n i <p UJ A x i) < Po J A x 2) 



(41) 



Using this representation, the two-particle density in the grand canonical ensemble is 
found to be 



R 2 {x,y) 



1 



_|_ e j3(u)i~ii) 



(42) 



Using similar manipulations as for the canonical partition function and the corre- 
sponding one-particle density it is possible to simplify the exact analytical expression for 
the connected two-particle correlation function in the canonical ensemble, 



N 



R 2 (x,y) = \MX)MV)\ 2 E (_l)'' e - 4 ^^+^-"(r+ S )+ 5 (r+,)(r+,+l)) JJ (1 _ <f). 

k,l r,s=0 p=N—r—s 



(43) 



This result can be used to compare the two ensembles but we will not address this question 
in this article. 
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4 Semiclassical Calculation 



In this section we calculate the microscopic spectral density and the two-point correlation 
function using semi-classical methods starting from the expressions for the grand canonical 
ensemble derived in previous section. We are thus interested in the region around x = 0. 
Because of the hard edge at x = we cannot simply do a WKB approximation by replacing 
the wave functions by plane waves but instead have to use Bessel functions. This follows 
immediately from the wave equation for x — > 0, 

y/xJ u (kx) = k 2 y[xJ v {kx). (44) 

Alternatively, one can exploit the asymptotic relation between Laguerre polynomials and 
Bessel functions. 



~d 2 x + 



4z/ 2 



Ax 2 



4.1 The Microscopic Spectral Density 

Taking into account the normalization of the Bessel functions, kdk\/ xx' J v (kx)J u (kx') = 
5(x — x'), to fix the constants in the integration measure, we arrive at the following ex- 
pression for the single particle density 

(p(x)) = fkxdk i + f^ uj2x2 . (45) 

The chemical potential is determined by the condition / dxp(x) = n. Since this integral is 
over all x, the use of the semi-classical expressions for the wave functions is not justified 
but instead we have to rely on the exact wave functions. In the limit, f3p 3> tv, the sum 
over % in fl3T|) can be replaced be an integral which can be performed analytically resulting 
in 

e (2u+2) P u-^ = 1 (46) 

In the limit n(3u ^> 1 the semi-classical expression for the spectral density is thus given 
by 

(P{X))=J Q kxdk - - e ^J^ X \n-2v-2)^ ■ ( 4? ) 



In this limit the semiclassical expression (|45|) leads to the correct value of the chemical 
potential. Below we will show for n(3uj ^> 1 and finite x (in units of the average level 
spacing) the term j3uj 2 x 2 can be neglected relative to An(3uj. 

An estimate for the average spectral density near zero but many level spacings away 
from x = is obtained by using the leading order asymptotic expansion of the Bessel 
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functions and calculating the integral in (|47| ) in the limit of a degenerate Fermi-gas. This 
results in 

P = (48) 

7T 



where k m y/Anu is the "radius" of the Fermi-sphere. At fixed h in the limit n > 1 we 
have 

uo « S 2 n and /3w w — . (49) 

2nn 

Using these results and invoking the Banks-Casher formula the parameter £ can be iden- 
tified as the chiral condensate, 

lim — = E. (50) 

We will now show that our approximations are self-consistent. The condition that 
we are close to the degenerate Fermi gas can be written as Sk/k <C 1. We thus have to 
impose the requirement that 

5k 1 h 

T~2^=4 <<L (51) 
The conditions kx ^> 1 and (3uj 2 x 2 <C Anf3u can be combined into 

— < ujx 2 < 4n, (52) 
4n 

or, in units of the average level spacing, x = u/p, the range of validity of the above 
asymptotic results is given by 

1 U 2 Tl 2 , s 

— < — — < An. (53) 
An An 

Because of the second inequality it is justified to neglect the term (3u 2 x 2 which we will 
do in the remainder of this section. 

By partial integration the expression (O) for the spectral density can be rewritten as 

= r ^^^(^-^(M^-i^^ (54) 

Jo A cosh f(/c 2 — Ancu) 

Using the Banks-Casher formula one finds that the chiral condensate depends on h. The 
leading order correction is given by 
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The corresponding spectral density will be denoted by p(h) = 2nS(/i)/7r. The microscopic 
spectral density is then given by 



Ps[U) 



V 1 ( u \ 

™2n£(/i) Pl 2n£(/i) j ' 



/l COsh 2 ^ r 1 - > 



E(/i) Jo 

where p° is the microscopic spectral density for the chGUE |53|, [54 



(56) 



P°s( u ) = 2 i J ^ M ) ~ 



(57) 



The interpretation is clear. The oscillations in the microscopic spectral density due to 
the Bessel functions are smeared out over a distance hu by the integration over k. The 
oscillations are thus visible up to a distance of u ~ 1/h. 



1.5 



7lp (71U) 



h = 0.3 (Exact) 

h = 0.3 (Approximate) 
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7TU 



Figure 1: The microscopic spectral density irp s (iru) versus nu for a value of h = 0.3. The 
exact result (full curve) and the approximate result (dashed curve) nearly coincide. 
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In the limit of small h the main contribution to the integral comes from the region 
around the Fermi-surface. In this limit we can derive an approximate formula correct to 
order h 2 at fixed uh. To this end we change integration variables in ( 56f) according to 
k = 1 + ht neglecting terms that are sub-leading in h. This results in 



Ps{u) 



cosh 2 (2t) 



(5f 



Next we Taylor expand the microscopic spectral density as follows 
p° s (u + hut) 



+ 
+ 

+ 
+ 



P» 

-hut 
2 

1 {hut) 2 

2 2! 

1 (hut) 3 

2 3! 

1 (hut) 4 

2 4! 



^ + Jv+lJu—! 



U V J V - X + — 

u 



A(J 2 V _ X - J, 



u+l. 



II 



u u u a 



(59) 



where the terms that will be neglected are denoted by /j, g^, hi, etc.. At small values 
of u these terms are of order h 2 u whereas for large u they are suppressed by order 1/u 2 
(notice the factor u in (|57D). By inspection one easily finds that the neglected terms are 
at most of order h 2 independent of the value of u. The leading order terms can be easily 
resummed to 

p° s (u + hut) = p° s (u) - ~J v (u)J u -x(u)[cos(2hut) - 1] + 0(h 2 ) + terms odd in t. (60) 



The integral over t in Q58) can be performed analytically resulting in 



Ps(u) = p° s (u) - -J v (u)J u -x(u) 



ix hu 



2 sinh(7r/ui/2" 



- 1 



+ 0(h 2 ). 



(61) 



In Fig. 1 we compare the exact expression for the microscopic spectral density (|56| ) to 
this approximate formula. We observe that even for a value of h as large as h = 0.3 the 
two results are very close. 



4.2 Two-Point Function 

In this subsection we derive a semi-classical expression for the two-point correlation func- 
tion. To this end the wave functions in the expression ( |4"2"| ) for the connected two-point 
correlation function are replaced by Bessel functions, 



R2(x,y) 



00 uju V^y J ^ kx "> J ^ ky ^ 

1-1- e /3k 2 -inf3uj 



(62) 
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By partial integration with respect to k the correlation function can be expressed as 

1 1 r x . /\ I err i >/,' // — ttk \ ' " 

R 2 (x,y) 



/> f°° u2 dk K(xnpk,ynpk) 



h Jo 



cosh'((F - l)/h) 



(63) 



where K(x, y) is the two-point kernel for the chGUE given by 

.xJ v+1 {x)J v {y) - yJ v {x)J v+1 {y) 



K(x,y) = y/xy- 



x 2 — y 2 



(64) 



To order h 2 the correlation function can be simplified to 



R 2 (x,y) 



P 



dt 



K(xirp(l + ht),yirp(l + ht) 



cosh^ 2t 



+ 0(h 2 ), 



(65) 



In the same way as for the one-point function we now will derive an approximate formula 
for the two-point function correct to order h 2 at fixed value of uh. This can be done 
conveniently by using the following summation formula 

Ju+i{x + xh)J l/ (y + yh) ± J v {x + xh) J v +i(y + yh) 

~ [J v+ x{x)J v {y) ± J u (x)J u+1 (y)] cos[(x + y)h] 

+ [J u (x)J u (y) T Ju+i(x)Ju+i(y)] sin[(x + y)h\. (66) 

This formula has been derived by means of a Taylor expansion and a subsequent resum- 
mation employing the following approximate derivative formulas 

dl n [J v+1 {kx)J u {ky) ± J v {kx)J„ +1 (ky)\ 

w (-l) n (x ± y) 2n [J v+l (kx)J v {ky) ± J v (kx)J u+1 (ky)l (67) 

and 



d^- x [J v+1 {kx)J v {ky) ± J v+1 {kx)J v {ky)\ 

w (-lT + \x±y) 2n - 1 [J u (x)J u (y)^J u+1 (x)J u+1 (y)]. (68) 

They have been obtained by means of recursion relations for Bessel functions neglecting 
terms that are suppressed by order 1/x or 1/y. One can easily show that the combined 
powers of x and y in the prefactor is always larger than the combined powers of x and 
y that have been neglected. Since only even terms in t contribute to the integral in (|65|) 
our final result for the two-point function, correct to order h 2 , is given by 



R 2 (x,y) 



n 2 hpy/xy 



+ 



% 2 hpy/xy 
8 



J v+ i(x7ip)J u (ynp) + J u (xirp)J u+1 (ynp) 
sinh((x + y)n 2 hp/4:) 

J v+1 (xirp)J u (ynp) - J u (xnp)J u+1 (yTip) 
sinh((x — y)ir 2 hp/4:) 



(69) 
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In the limit x, y ^> \x — y\ the analytical result for the two-point function of the model 
of Moshe, Neuberger and Shapiro [40] is recovered from the leading order asymptotic 
expansion of the Bessel functions. For unitary invariant ensembles, it can be shown 
that the result of |^0| for critical statistics and Wigner-Dyson statistics are the only two 
possibilities | 55| . At this moment it is not clear whether this argument can be extended 
to the chiral unitary ensembles as well. 




Figure 2: The number variance S 2 (L) versus L for a value of h = 0.3. The exact result 
(full curve) and the approximate result (dashed curve) nearly coincide. 



The number variance of the eigenvalues near x = is obtained by integrating the 
two-point correlation function including the self-correlations 

r L/p(h) f L/p{h) 

S 2 (L)=/ dx dy[5(x-y)(p(x)} + R 2 (x,y)}. (70) 

Jo Jo 

We study its asymptotic behavior in the limit L — > oo. Starting from the expressions 
(U) and (|62|) , S 2 (L) can be simplified by means of an orthogonality relation for Bessel 
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functions. To leading order in h we find 

S 2 (L) = / xdx 

Jo Jo 



L /PW roc Jl( kx \ 

' kdk- 



4 cosh 2 | (k 2 - 4nuj) 



h L. (71) 



p2%pk 4 

The same asymptotic result can be derived from the analytical result ( |B"9"D (In this case 
the average spectral density does not depend on h (see |H).). Such linear term, first 
proposed in ||31|| , is believed to be characteristic for universal critical statistics [|37| valid 



at the mobility edge and has been related to the multifractality of the wave functions 



5j||5^£§|59j. 

In Fig. 2 we compare the number variance derived from the approximate analytical 
result ( |69|) and from the exact result (|62|) . Clearly, even for a value of h as large as 0.3 
the two curve are barely distinguishable. 

The asymptotic linear behavior of the number variance seems to be contradicted by 
the sum rule 

roo 

dy[5(x - y)(p(x)) + R 2 (x, y)\ = 0. (72) 



The resolution of this paradox |4T], [6(| is probably best illustrated by considering the 
Poisson ensemble for n uncorrelated eigenvalues with average spacing p. To satisfy the 
sum-rule we have that R^{x, y) = p/n instead of zero for uncorrelated eigenvalues resulting 
in the number variance ^(L) = L—L 2 /n. We conclude that an asymptotic linear behavior 
is possible if the thermodynamic limit is taken before the limit L — > oo. 

To make contact with the partially quenched effective partition function, for which 
the number of subsequent eigenvalues around zero that are correlated according to the 
chGUE scales as y/n @, |5|, [], 0, we have to scale h as h — > h/ y/n. In this limit microscopic 
universality Jl|, [8|, |9], [10], 0, [12|] is recovered for the interpolating chiral unitary ensemble. 



5 Comparison with Instanton Simulations 

Spectral correlations have been studied in great detail for both lattice QCD simulations 
0, |T^, 0, [19], [T7|, |21| and instanton-liquid simulations (see for a recent review). 
In lattice QCD they were studied by means of the disconnected scalar susceptibility |22|j , 
and complete agreement with partially quenched chiral perturbation theory was found. In 
particular, it was shown that the number of subsequent eigenvalues around zero described 
by chRMT scales as F 2 \/V . In instanton simulations a weaker volume dependence of 
the number of such eigenvalues was observed suggesting an approach to a critical point 
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o N=128 




L 

Figure 3: The number variance S 2 (L) versus L for Dirac eigenvalues of instanton gauge 
field configurations (points) and our interpolating random matrix model. The total num- 
ber of instantons and the parameter h of our model are given in the label of the figure. 

similar to a localization transition. Indeed the multifractality index of the fermionic wave 
functions was found to be nonzero. We thus compare the instanton data with the model 
in previous section at fixed value of the parameter h. Results for the number variance, 
S 2 (L) versus L are shown in Fig. 3. The closed and open circles represent results |J for 
the eigenvalues of the Dirac operator with field configurations given by an ensemble of 
instantons and an equal number of anti- instantons with a total density of 1 /m 4 . The 
total number of (anti-)instantons is given in the label of the figure. In the same figure we 
show the result for the chGUE (dotted curve) and results for the model (^3) for h = 0.18 
(full curve) and h = 0.23 (dashed curve). We observe that both the slope and the range 
of agreement with the chGUE curve only shows a week volume dependence. Outside this 
domain the data show a linear L-dependence. Both features are nicely reproduced by 
the critical Random Matrix Model. The slightly positive curvature of the instanton data 
might be a remnant of the L 2 dependence predicted for the Altshuler-Shklovsky domain 
0. 
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6 Conclusions 



We have analyzed a chiral random matrix model that interpolates between the chGUE and 
the chiral Poisson ensemble. This model is a generalization of a model originally proposed 
by Moshe, Neuberger and Shapiro |^D[. It has been mapped onto a gas of non- interacting 
fermions and was solved by means of statistical mechanics methods. To leading order 
in the deviation from the chGUE we have obtained compact analytical expressions for 
the microscopic spectral density and the two-point level correlation function. We have 
shown that this critical chiral random matrix model provides a good description of the 
level correlations of Dirac eigenvalues for gauge field configurations given by a liquid of 
instantons. 

The number variance of the critical chiral random matrix model shows a linear in- 
dependence for a large L (in units of the average level spacing) whereas it coincides with 
the chGUE result for small values of L. The characteristic feature is that the transition 
point between these two domains is stable in the thermodynamic limit. This situation 
is very different for a non-linear cr-model description of disordered systems where this 
transition point or the Thouless energy is determined by the competition between the 
mass term and the kinetic term. In that case one finds the scaling behavior E c ~ D / ' L 2 S 
with D the diffusion constant and L s the linear size of the sample. The theoretical reason 
for a scale independent dimensionless conductance (i.e. the Thouless energy in units of 
the average level spacing) is that the localization length diverges at a critical value of 
the disorder. In the approach to this limit the diffusion constant has to become scale 
dependent. If E c , in units of the average level spacing, becomes scale independent the 
diffusion constant has to be scale dependent leading to a multi-fractal scaling of the wave 
functions. 

The weak volume dependence and the linear number variance observed in correlations 
of eigenvalues of the QCD Dirac operator with instanton liquid gauge field configurations 
suggests that we are dealing with a critical system close to a localization transition. 
Indeed the same numerical simulations suggest a small nonzero multifractality index of 
the wave-functions. On the other hand, the dimensionless conductance found in lattice 
QCD simulations scales according to our expectations from chiral perturbation theory. At 
this moment we do not have a good explanation for this discrepancy. It could simply be 
that the expected scaling behavior is only recovered for much large volumes in instanton 
simulations. Indeed, a very slow approach to the thermodynamic limit has been found for 
other quantities such a quenched chiral logarithms. Clearly, more work has to be done to 
resolve this issue. 
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